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We present a nonperturbative renormalization-group approach to the Bose-Hubbard model. By 
taking as initial condition of the renormalization-group flow the (local) limit of decoupled sites, 
we take into account both local and long-distance fluctuations in a nontrivial way. This approach 
yields a phase diagram in very good quantitative agreement with quantum Monte Carlo simulations, 
and reproduces the two universality classes of the superfluid-Mott-insulator transition. The critical 
behavior near the multicritical points, where the transition takes place at constant density, agrees 
with the original predictions of Fisher et al. [Phys. Rev. B 40, 546 (1989)] based on simple scaling 
arguments. At a generic transition point, the critical behavior is mean-field like with logarithmic 
corrections in two dimensions. In the weakly-correlated superfiuid phase (far away from the Mott 
insulating phase), the renormalization-group flow is controlled by the Bogoliubov fixed point down 
to a characteristic (Ginzburg) momentum scale fee which is much smaller than the inverse healing 
length kh- In the vicinity of the multicritical points, when the density is commensurate, we identify 
a sharp crossover from a weakly- to a strongly-correlated superfiuid phase where the condensate 
density and the superfiuid stiffness are strongly suppressed and both ka and kh are of the order of 
the inverse lattice spacing. 
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I. INTRODUCTION 

In the last two decades, the nonperturbative renormal- 
ization group (NPRG) approach has been successfully 
applied to many areas of physics, 1,2 from high-energy 
physics to statistical and condensed-matter physics. It 
has proven to be a powerful tool to study not only 
the low-energy long-distance properties in the vicinity 
of second-order phase transitions but also non-universal 
quantities. In particular, the NPRG approach has been 
implemented in lattice models and used to compute the 
transition temperature and the magnetization in classi- 
cal spin models (Ising, XY and Hcisenberg models). ' This 
implementation of the NPRG is referred to as the lattice 
NPRG. 

The strategy of the NPRG is to build a family of 
models indexed by a momentum scale k, such that fluc- 
tuations are smoothly taken into account as k is low- 
ered from a microscopic scale A down to 0. In practice 
this is achieved by adding to the action S of the sys- 
tem an infrared regulator term AS}, which vanishes for 
k = 0. For a scalar field theory, the regulator term is a 
mass-like term ASk[f] = 2 ¥ , -q-^fc( { l)¥'q> where the 
cutoff function i2fc(q) is chosen such that iifc(q) ~ fc 2 
f° r IqI ^ & an d -Rfc(q) ~ for |q| > fc, which effec- 
tively suppresses the low-energy modes |q| < fc. One can 
then define a scale-dependent partition function Zk[J] 
and a scale-dependent effective action rj-[^>] defined as 
a slightly modified Legendre transform (see Sec. II A for 
the precise definition) of — In Z^ [J] . Here J is an ex- 
ternal source which couples linearly to the (p field and 
4>(t) = 5 In Zk[J]/SJ(v). In the standard implementation 
of the NPRG, at the microscopic scale fc = A, all fluctu- 
ations are frozen by the ASa term so that Ta[4>] = S[4>] 



as in Landau's (mean-field) theory of phase transitions. 
The effective action of the original model is obtained for 
fc = (ASk=o — 0) and can be determined by (approxi- 
mately) solving the RG equation satisfied by T^. ' 

The lattice NPRG differs from the standard implemen- 
tation in the initial condition. The cutoff function (q) 
is chosen such that at the microscopic scale fc = A the 
action S + AS a corresponds to the local limit of decou- 
pled sites. Local fluctuations are therefore included from 
the very beginning of the RG procedure. The intersite 
coupling is then gradually restored as fc decreases from A 
down to 0. In the low-energy limit fc -C A, i?fc(q) acts as 
an infrared regulator suppressing fluctuations with mo- 
menta |q| < fc. The lattice NPRG is then equivalent to 
the standard NPRG and yields identical results for the 
critical properties. The hallmark of the lattice NPRG is 
thus to take into account both local and critical fluctua- 
tions in a nontrivial way. 

In this paper, we present a NPRG study of the Bose- 
Hubbard model at zero temperature and in dimension 
d = 2 or d = 3. This model has been intensively stud- 
ied in the last years following the experimental obser- 
vation of the superfluid-Mott-insulator transition of an 
ultracold bosonic gas in an optical lattice. '~ lS Phase dia- 
gram and thermodynamic quantities are known from the 
numerically exact lattice quantum Monte Carlo (QMC) 
simulations."' On the other hand few studies have ad- 
dressed the critical behavior at the superfluid-Mott- 
insulator transition, and most of our understanding 
goes back to the seminal work of Fisher et al. 4 

The standard NPRG scheme does not capture 
the superfluid-Mott-insulator transition in the Bose- 
Hubbard model. The reason is that near the transition 
the mean-field solution is too far away from the actual 
state of the system to provide a reliable initial condi- 
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tion for the NPRG procedure. The two-pole structure of 
the local (on-site) single-particle propagator is crucial for 
the very existence of the transition (see, e.g., Ref. 4). It 
is however impossible to reproduce this structure from 
a RG approach starting from the mean-field (Bogoli- 
ubov) theory within standard approximations of the RG 
equation satisfied by the effective action This pre- 
vents a straightforward generalization of recent NPRG 
studies 12-18 of interacting bosons to the Bose-Hubbard 
model. 

By contrast the lattice NPRG, which takes into 
account local fluctuations, is able to describe the 
superfluid-Mott-insulator transition. 1 ' 1 Since the start- 
ing action S + A.S\ is purely local, this approach is to 
some extent reminiscent of various t/U expansions of the 
Bose-Hubbard model. 20-2 ' Moreover, the lattice NPRG 
is not restricted to the computation of thermodynamic 
quantities and allows us to study the critical behavior 
at the superfluid-Mott-insulator transition and compare 
with the predictions of Fisher et al. 1 based on scaling 
arguments. 

In addition to the phase diagram and the critical be- 
havior at the superfluid-Mott-insulator transition, the 
NPRG approach can also address the superfluid phase. 
Deep in the superfluid phase, localization effects are neg- 
ligible and we expect the Bogoliubov theory to provide 
a good description of the system. However, even in 
this weak correlation limit, it is known that the Bogoli- 
ubov approximation breaks down below a characteristic 
(Ginzburg) momentum scale kc- In perturbation the- 
ory about the Bogoliubov approximation, the Ginzburg 
scale manifests itself by the appearance of infrared diver- 
gences below three dimensions (d < 3). Although these 
divergences cancel out in local gauge invariant quanti- 
ties (condensate density, sound mode velocity, etc.), 28 ~' !1 
they do have a physical origin: they result from the cou- 
pling between longitudinal and transverse (phase) fluc- 
tuations and reflect the divergence of the longitudinal 
susceptibility ' - a general phenomenon in systems 
with a continuous broken symmetry. The normal and 
anomalous self-energies, S n (q, uS) and £ an (q, w), are non- 
analytic functions of q and uj when |q|, \uj\/c <C &g an d 
d < 3 (c denotes the velocity of the sound mode), while 
£ an (0, 0) vanishes, in marked contrast with the Bogoli- 
ubov approximation where the linear spectrum and the 
superfluidity rely on a finite value of the anomalous self- 
energy. A weakly-correlated superfluid is defined by the 
condition kc *C kh where the healing scale kh is the in- 
verse of the healing length £/, = fe- 1 . 14 ' 39 . 41 ! 42 I n this 
case, the Bogoliubov theory applies to a large part of 
the spectrum where the dispersion is linear (|q| <C kh) 
and breaks down only at very low momenta |q| <C kc- 
The Goldstone regime |q| <C k G , dominated by phase 
fluctuations, is conveniently described by Popov's hydro- 
dynamic theory (free of infrared divergences) based on 
a density-phase representation of the boson field tp — 
yjne 1 ® . 43-45 The NPRG approach yields a unified descrip- 
tion of superfluidity which includes both Bogoliubov the- 



ory (valid for |q| ^> kc) and Popov's hydrodynamic ap- 
proach (valid for |q| < k h ). 13 ' 14 ' 39 The Bose-Hubbard 
model gives us the opportunity to understand the fate of 
the weakly-correlated superfluid phase as we increase the 
strength of the interactions and move closer to the Mott 
insulating phase in the phase diagram. 

The paper is organized as follows. In Sec. II, we de- 
rive the lattice NPRG formalism for the Bose-Hubbard 
model. We introduce the scale-dependent effective ac- 
tion Tk and compute its initial value T\ by solving the 
single-site Bose-Hubbard model. We show that T\ re- 
produces the result of the strong-coupling random-phase 
approximation (RPA). 46-50 We also discuss the approxi- 
mations used to solve the flow equation satisfied by r^. 
The phase diagram obtained from the NPRG equations 
is in very good quantitative agreement with the QMC 
results (Sec. III). Furthermore, the critical behavior de- 
rived from the NPRG analysis is in complete agreement 
with the predictions of Fisher et al. based on scaling ar- 
guments (Sec. IV). 4 We find multicritical points in the 
universality class of the (d + l)-dimensional XY model 
where the transition takes place at constant density. The 
XY critical behavior is observed in the Mott gap, the 
condensate density, the compressibility and the super- 
fluid stiffness when a multicritical point is approached 
at constant chemical potential by varying the ratio t/U 
between the hopping amplitude and the local repulsion 
between particles. At a generic transition point, we ob- 
serve mean-field behavior, with logarithmic corrections 
in dimension d = 2 (corresponding to the upper critical 
dimension). The superfluid phase is discussed in Sec. V. 
In the dilute limit, the renormalization-group flow is con- 
trolled by the Bogoliubov fixed point down to a charac- 
teristic (Ginzburg) momentum scale fee which is much 
smaller than the inverse healing length kh- The Gold- 
stone regime k <§; kc, dominated by phase fluctuations, 
is characterized by a (relativistic) Lorcntz invariance of 
the effective action T^. ■ In the vicinity of the mul- 
ticritical points, when the density is commensurate, we 
identify a sharp crossover from a weakly- to a strongly- 
correlated superfluid phase where the condensate density 
and the superfluid stiffness are strongly suppressed and 
both kc and kh are of the order of the inverse lattice 
spacing. The main results are summarized in Sec. VI. 

II. LATTICE NPRG 

The Bose-Hubbard model on a c?-dimensional hypercu- 
bic lattice is defined by the (Euclidean) action 

5 = J Q dT { E " M)Vr + ^(VV>r) 2_ 

-t£(l/S^+c.c.)}, (1) 

where ip r (r) is a complex field and r £ [0, /3] an imaginary 
time with f3 — > oo the inverse temperature, {r} denotes 
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FIG. 1. (Color online) Effective (bare) dispersion t q + Rk(q) 
for k = A, < k < A and k = with cutoff function (4). The 
(green) dashed line shows the bare dispersion t q = — 2icosq. 
d — 1 and t is taken as the energy unit. 



the N sites of the lattice. U is the on-site repulsion, t the 
hopping amplitude between nearest-neighbor sites (r, r') 
and fj, the chemical potential. In the following we will 
sometimes write the boson field 



i'r = — S (^lr + #2r) 



(2) 



in terms of two real fields ipi r and tp^r- 

We set h — ks = 1 and take the lattice spacing as the 
unit length throughout the paper. 



A. Scale-dependent effective action 

Following the general strategy of the NPRG, we con- 
sider a family of models with action S k = S + AS k in- 
dexed by a momentum scale k varying from a microscopic 
scale A down to 0. The regulator term is defined by 



AS k 



Gh-2^V>q-Rfc(q)V>q 



(3) 



where ipq is the Fourier transform of ip r and the sum 
over q runs over the first Brillouin zone ] — 7T, ir] d of the 
reciprocal lattice. The cutoff function -Rfc(q) modifies 
the bare dispersion i q = — 2t Yli=i cos 1i °f tne bosons. 
i?A(q) is chosen such that the effective (bare) dispersion 



The scale-dependent effective action 

T k [<f ,<t>]= -lnZ k [J*, J] + / drV(j; 

Jo 



+ C.C.) 



(5) 



is defined as a (slightly modified) Legendre transform 
which includes the explicit subtraction of AS k [4>* , (j)]. 
Here Z k [J* , J] is the partition function obtained from 
the action S + ASk, J r a- complex external source which 
couples linearly to the bosonic field ip r , and 

S\nZ k [J*,J] S\nZ k [J*,J] 
<M T ) = — c tw-\ — ' ^r( T ) = rrn — (6) 



£J r *(r) 



<JJr(r) 



the superfluid order parameter. The variation of 
the effective action with k is governed by Wetterich's 
equation, 01 



0k W 



1 



T^d k R k (r fe 



(2) 



(7) 



(2) 

where IT is the second-order functional derivative of 
Tfc. In Fourier space, the trace in (7) involves a sum over 
momenta and frequencies as well as the two components 
of the complex field 4>. 

We are primarily interested in two quantities. The first 
one is the effective potential defined by 



V k {n) 



1 



(8) 



const 



where <f> is a constant (uniform and time-independent) 
field. The U(l) symmetry of the action implies that V k (n) 
is a function of n = \tf)\ 2 . Its minimum determines the 
condensate density no tk and the thermodynamic poten- 
tial (per site) Vo,fe = Vfc(no,fc) in the equilibrium state. 
The second quantity of interest is the two-point vertex 



.(2) 
k,ij 



(r — r', r — t' 



(9) 



tq + R\(q) vanishes/' The action S A = S + AS a then w hi c h determines the one-particle propagator G k 
corresponds to the local limit of decoupled sites (vanish- _-p( 2 ) -1 
ing hopping amplitude). By choosing Ra{<1) + tq = 0, 
rather than Ra(<1) +t q = 2dt as in Ref. 3, we ensure that 
R\(q_) does not modify the chemical potential but only 
the kinetic energy. 

In practice, we choose the cutoff function 

R k (d) = -Z 4 , fc e fc sgn(< q )(l - y q )0(l - y q ), (4) 



with A = V2d, e k = tk 2 , y q = (2dt - |t q |)/cfe and Q(x) 
the step function (see Fig. 1). The fc-dependent constant 
Za,r is defined below (Za,a — !)■ Since R k= o{<l) = 0, the 
action 5^=0 coincides with the action (1) of the original 
model. For small k, the function R k (q) gives a mass ~ k 2 
to the low-energy modes |q| < k and acts as an infrared 
regulator as in the standard NPRG scheme. 1 ' 2 



Here the indices i, j refer to the real and imagi- 
nary parts of 4> [see Eq. (2)]. Because of the U(l) symme- 
try of the action (1), the two-point vertex in a constant 
field takes the form 14 

r ifi-fe0) = S i jT A , k (q;n)+(t) i (j} j r Btk (q;n)+e l: jT c ^(q;n) 

(10) 

in Fourier space, where q = (q, iu>), uj is a Matsubara 
frequency and the antisymmetric tensor. For q = 

(2) 

0, we can relate 1^ to the derivative of the effective 
potential, 



ri 2 i(9 = o ;i 



(11) 
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so that 

T Ark (q = Q;n) = Vl(n), 

T Btk (q = 0;n) = V k "(n), (12) 

r Cjfc (g = 0;n) =0. 

Parity and time-reversal invariance imply 

T A . k (q;n) = T Atk (-q;n) = T Aik (q, -iu; n), 

r_B,fe(g; n) = T B>k (-q; n) = r B!fe (q, -iuj; n), (13) 

Fc,k(q; n) = -r C)fe (-g; n) = -r C)fe (q, -iu; n). 

(2) — 1 

The one-particle propagator G k = — Fj. ; can be 
written in a form analogous to (10) or in terms of its 
longitudinal and transverse components, 



Gk,ij(q; <t>) 



G k ,n(q;n) + 5ij 



2n 

+ £ijG k it(q;n), 



) G kttt (q;n) 



(14) 



where 



G k ,n(q;n) 



T Ajk (q;n) 
D k (q;n) 



n i x F A _ k (q;n) + 2nT Byk (q;n) 
G kM (q;n) = : n yz7Z\ ' ' ( 15 ) 



D k (q;n) 



G k ,\t{q\n) 



D k (q;n) ' 



with D k = T\ k + 2nT Aik T Byk + T 2 ck . Note that the 
single-particle propagator entering the flow equation (7) 

is defined by — (1^ + -Rfc) -1 , which is the propagator 
associated with the true Legendre transform, rather than 



(2)-l 



-r 



B. Initial conditions 

Since the action S + AS\ = S\ oc corresponds to the 
local limit, the initial value of the effective action reads 

r A [0*,0] = r loc / drV^( q )v(q), (16) 

J ° q 

where 



loc 



r0 „ 

4] = -lnZio C [J*,J]+ / rfrV(j;0 r + c.c.) 
Jo 

(17) 

is the Legendre transform of the thermodynamic poten- 
tial — In Z\ oc [J* , J] in the local limit. In Eq. (17), J is re- 
lated to (j> by the relation r (r) = 5\nZ\ oc {J* , J]/5J*(t) 
and Z\ oc is the partition function obtained from S\ oc . 

It is not possible to compute the functional Ii oc \(jf , <fi] 
for arbitrary time-dependent fields. One can however 
easily obtain the effective potential Vi oc (n) and the two- 

(2) 

point vertex rj in a time-independent field (j). These 




FIG. 2. (Color online) (Top) Superfluid order parameter <f> 
vs external source J (here assumed real) in the local limit for 
various values of the chemical potential fj,. (Bottom) Effective 
potential Vj c("-)- 



quantities are sufficient to specify the initial conditions 
of the flow within the approximations that we consider 
in Sec. HE. 

To obtain V\ oc and T\ oc in a time-independent field, it is 
sufficient to consider a single site with time-independent 
complex external source J. The corresponding Hamilto- 
nian reads 

H=-p,n+~h(h-l)- J*b- jtf, (18) 

where w (b) is a creation (annihilation) operator and n = 
b*b. In the basis {|?n)} [n|m) = m\m) with m integer], 
the Hamiltonian is represented by a tridiagonal matrix, 

(m\H\rn!) = & m ,m' —fJ.ni -\ m(m — 1) 



i J* \Jm + 1 — S„ 



(19) 

which can be numerically diagonalized in the truncated 
Hilbert space m < m max . The low-energy eigenstates 
are independent of m max if the latter is large enough. If 
we denote by {\a),E a } the source-dependent eigenstates 
and eigenvalues - with {|0),i?o} the ground state - we 
obtain the superfluid order parameter 

BE, 8E (2Q) 



d,J* ' 

and the effective potential 



8J 



Vi oc (n) = E Q + J*4> + J(f> 



(21) 
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FIG. 3. (Color online) [T loCzA (iuj; n) - V{ oc {n)]/U vs u/U for 
various values of n. fj, = 0.2U (left) and fx = (V2 — 1)17 
(right). The dotted lines show the large frequency limit (25). 
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FIG. 4. (Color online) rioc,c(*w; n)/U vs w/f7 for vari- 
ous values of n. y. = -0.2(7, 0.2(7, (V§ - 1)(7 and 0.6(7 
(from top left to bottom right). In the large frequency limit, 
T\ oc ,c(iu n ; n) = uj [Eqs. (25)]. 



(n = \tfi\ 2 ) in the zero-temperature limit j3 — > oo. Figure 2 
shows the superfluid order parameter <p as a function 
of the external source J, and the local effective poten- 
tial V\ oc (n) obtained by numerically inverting (20). The 
special case where the ground state in the local limit is 
degenerate for J = (fi/U integer) is discussed in Ap- 
pendix A. 

(21 

To determine the two-point vertex Ty oc , we start 
from the (source-dependent) normal and anomalous local 
Green functions 



Hamiltonian, 



\(0\b\a)\< 



iuj + E a — Eq iuj + Eq — E a 



From the relation T^ 2 ) = — G 1 , we obtain 



1 



(23) 



Ti OCtA (iuj;n) = - — [G n (iw) + G n (-iw) + 2G an (iw)], 



rioc,s(«o;;n) = 



Gan(^) 

nD ' 



r loC)C («w;n) = —[G n (ioj) - G n (-«*;)], 



(24) 



where D — G n (iu})G n (—iuj) — G an (iu;) 2 . T^ 2 ^ is expressed 
in terms of the condensate density n (rather than the 
external source J) by inverting (20). 

The large frequency limit of the two-point vertex is 
given by 52 



lim Ti OCiA {iu; n) = -fi- U{ip 2 ) + 2U (^*ip r ), 
Ujj 2 ) 

lim r locB (iu;;n) = U — — , 
|w|->oo n 

lim r loC; c(iw; n) = u, 

\uj\—>oo 



(25) 



where (ip r ) = 4> and (ipr) are assumed real (which corre- 
sponds to a real external source J), 

-G n (r = O-) + |0| 2 

dui „ , . , iuj0 + 



2tt 



G n {iuj)e 1 ' 



\<t>Y 



(26) 



and 



{€) =-G an (r = O-) + 2 



2tt 



G an (iw)e l 



(27) 



The asymptotic forms (25) are reached for |w| 3> U. 
rioc,y4(*w;n) and ri oc c^iu;; n) are shown in Figs. 3 and 
4. 



1. Large-field limit 



G n {T) = -(T T b{r)Pm + \(b)\\ 
G an (r) = -(T r S(r)S(0)) + (S) 2 , 



(22) 



where W\t) — e rH b^ "e~ rH and T T is a time-ordering 
operator. The Fourier transforms G n (iu>) and G an (iw) 
are easily expressed in terms of the eigenstates \a) of the 



To obtain the large-field limit of the effective poten- 
tial V\ oc (n), we must compute the local partition func- 
tion Z\ oc [J* , J] for I J I — > 00. Using a loop expansion 
about the saddle-point approximation (see Appendix B), 
we find 

VIocM = -fin + ^n 2 + O{n ), (28) 



where fi = fi + U(l- y/3/2). 
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2. Strong- coupling RPA 

The initial effective action Ta [Eq. (16)] treats the local 
fluctuations exactly but includes the intersite hopping 
term at the mean-field level, thus reproducing the strong- 
coupling RPA. 46-50 The effective potential reads 

V A (n) = V loc (n) - 2dtn, (29) 

while the two-point vertex takes the RPA-like form 

r$j (<?; n) = r™ ti . n) + s id t^. (so) 

Expanding Va(ji) about n — 0, we find 

V A (n) = ^(Oj+Jrg^^iw = 0;n = 0) - 2eft] n+G(n 2 ), 

(31) 

where 

rSc,«(*w = 0; n = 0) = -G a (iu> = 0; n = 0)" 1 (32) 
is determined by the local Green function 

G a {iu;n = = - — ■ — — ■ — — 

ilo + fi — U ri\ oc vjj + u — U (nio C — 1) 

(33) 

for vanishing source (J* = J = 0). Here ni oc is the 
number of bosons per site in the local limit: 7ii oc — 1 < 
fi/U < h\ oc if (J, > and ri ioc = if /i < 0. The ground 
state is a Mott insulator as long as V^(0) > 0. Thus the 
transition to the superfluid state is determined by the 
criterion V A (0) — 0, i.e. 



G n (iu = 0; n = + 2dt = 0, 



(34) 



which reproduces the mean-field (or strong-coupling 
RPA) phase diagram. 1 Equation (34) can also be ob- 

uj = 0;n = 0), 



(2) 

tained from the condition detT A (q 



which signals the appearance of a pole at zero mo- 
mentum and frequency in the one-particle propagator 

G A = -rf-\ 

In the strong-coupling RPA, the condensate density no 
in the superfluid phase is determined by 



Va(«o) = V£ c (no) - 2dt = 0. 



(35) 



The hopping amplitude t acts as a source term for the 
local potential Vi oc {n). For t/U 3> 1, i.e. deep in the 
superfluid phase, the source term is large and we are ef- 
fectively in the large field limit discussed in Sec. II Bl. 
From Eqs. (16), and the fact that p, + 2dt ~ u + 2dt 
when t/U ^> 1, we then obtain V\{n) ~ -^S[(j)* , <f>] (with 

1 0| 2 = n), which is nothing but the result of the Bogoli- 
ubov approximation. 11 ' ' ' The strong-coupling RPA re- 
duces to the Bogoliubov theory in the limit t/U » 1. 
Table I compares the initial conditions given by the Bo- 
goliubov theory and the strong-coupling RPA. 

It should be noted that the true Legendre transform 
is ri oc for k = A since the action S\ = S + ASa = Si oc 



Bogoliubov 



RPA 



Va(ti) -(fj, + 2dt)n+ fn 2 Vi oc (n) - 2dtn 
FA,A(q;n) —fi + Un + tq T loCi A{iLo;n) + t q 
TB,A(q;n) U T loCt B(ioj;n) 



TABLE I. Initial conditions given by the Bogoliubov approx- 
imation and the strong-coupling RPA. In the weak-coupling 
limit the two approximations become equivalent. 



is local. The lattice NPRG is an expansion about the 
local limit. The scale-dependent effective action Tfe is 
however the right quantity to consider to analyze the 
physical properties of the system with action S (with- 
out the regulator term). In T^, the regulator term ASk 
is compensated, in a mean-field manner, by subtracting 
ASk [</>*, <j>] from the true Legendre transform. It follows 
that the physical quantities at scale k, such as the con- 
densate density n k , are obtained from rather than 
from the true Legendre transform. 



C. Gauge invariance and Ward identities 

The invariance of the action S + AS k in the local 
(time-dependent) gauge transformation ip r — > e la ip T , 
tp* — > e~ za ipr and /x — > /x + id T a imposes important 
constraints on the effective action T^. In the superfluid 
phase, this implies that the two-point vertex satisfies the 
Ward identities 



_d_ 
Duj 2 



d 2 v k 



dndfi 
1 



q=0 



2n , k dfi 2 



(36) 



where the effective potential Vk(n, fi) is considered as a 
function of both n and u, the condensate density n^_ k = 
n o,k(l^) being then defined by 



dV k (n,fi) 
dn 



= 0. 



(37) 



Since Eq. (37) is valid for any /x, we deduce the relation 

(38) 



d dV k 




d 2 v k 


d 2 V k 


dn 0}k 


dji dn 


no, is 


dndfi 


«0,fc 


dn 2 


n , k d » 







which will be used below together with (36). 

In the Mott insulator (no,k = 0), the Ward identi- 
ties (36) become 



d_ 

duj 



^c,k{q\n = 0) 



d 2 V k 



d 2 V , k d 2 V ok 



q=0 
2l 



dy 2 



dndfi 
= 0, 



n=0 



(39) 
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(40) 



which implies that the compressibility 54 

_ dfik __ rf 2 Vb,fc 
dp dp 2 

vanishes, n k = —dVo k/dp denotes the boson density 
(mean boson number per site). 

D. Derivative expansion and infrared behavior 

The low-energy behavior of the system is best under- 
stood from a derivative expansion of the two-point ver- 
tex. Since the cutoff function (4) acts as an infrared 

(2) 

regulator, 1^ (q; n) is a regular function of q for q — > 0. 
In the infrared limit, we can therefore use the derivative 
expansion 



^A.k{q;n) 

Fc,k(q;n) 



ZA,k{n)tq' 
Zc,k{n)ui, 



+ V A , k (n)u J 2 + Vl(n), 



(41) 



{V k (n) = dVk/dn, etc.) in agreement with the symme- 
try properties (13). For the following discussion, it is 
convenient to introduce 



S k = 



dVk 
dn 



Xk = 



d 2 V k 
dn 2 



(42) 



with 5). vanishing in the superfluid phase. If the spectrum 
is gapped, Eqs. (41) will always be valid for energy scales 
below the gap. Otherwise their validity requires |q| < A; 
and \uj\ < cu k where tu k is the lowest excitation energy 
for |q| ~ k (see Sec. HE). 



1. Superfluid phase 

Using (37) and (38), we can rewrite the Ward identi- 
ties (36) as 



Zc,k(no,k) — X k 
VA,k{n>o,k) = 



dn 0yk 
dp ' 
1 d 2 V k 



2n 0jfc dp 2 



(43) 



no,k 



while the compressibility (40) is expressed as 



d 2 V ( 



0,k 



dp 2 



d 2 V . k 



dndp 



«0,fc 



dn 0}k 
dp 



2n ok V A h + 



7 2 
^Ck 

At 



(44) 



The superfluid stiffness p Syk , defined as the rigidity 
with respect to a twist of the phase of the order param- 
eter, can be obtained from the transverse part of the 
two-point vertex, 14 



r A , fc (q, W = 0;n . fc ) = ^q 2 
2n k 



0), 



(45) 



which leads to 



Ps,k = 2tZ A , k (n , k )n , k - 



(46) 



The excitation spectrum is given by the zeros of the 

(2) 

determinant of the 2x2 matrix T k (q;no ik ) (after ana- 
lytical continuation icu — > cu + i0 + ), 



detrf( 9 ) 



^A,k(q) l^A,k(q) + 2n 0ife r Bifc (q)] + T c ^ k (q) 2 
~ 2X k n 0tk (Z A . k t(i 2 + V A . k u 2 ) + {Zc,ku) 2 

(47) 

(all quantities are evaluated for n — n o k ) for |q|, \uj\ — ¥ 0. 
This equation yields a gapless (Goldstone) mode ui = 
Cfc|q| with a velocity 



Cfe 



ZA,k(no,k)t 



VA,k( n 0,k) 
\ 1/2 

Ps,k 
K k 



+ Z c , k {n a ^ k ) 2 / (2A fc n 0jfc ) 



1/2 



(48) 



which can be expressed in terms of the compressibil- 
ity and superfluid stiffness. The existence of a gap- 
less mode is a consequence of the Hugenholtz-Pines 
theorem which, in our formalism, reads ' 



^A,k(q = 0;^o,fe) = v k (no,k) 



= o. 



(49) 



Equations (46,48) are identical to those obtained in con- 
tinuum models if we identify l/2t with the (effective) 
mass to of the bosons in the lattice potential. From 
detr^ 2) (q,w + z0+) =0, we also obtain a gapped mode, 
with a gap which is however larger than lu^ ~ c k k, and 
therefore beyond the domain of validity of the derivative 
expansion (|q|, |w|/cfe <C k). The existence of two modes 

( 2") 

in the superfluid phase follows from det ' (q) being of 
order w 4 . Pushing the derivative expansion to higher or- 
der in u! 2 would yield additional modes. These modes are 
not in the domain of validity of the derivative expansion 
and do not show up in the spectral function. 14 ' 18 



2. Mott insulator 

In the Mott insulator (no,fc = 0), the Ward identi- 
ties (39) yield 



^C.k 



dS k 
dp 



(50) 



Since G alltk (q) = 0, the excitation spectrum is obtained 
from G~ j.(q) = -T A , k (q) + iT c ,k{q) = after analytical 
continuation iu — > u + i0 + . This gives two gapped modes 



w±(q) 



Zc,k _j_ 1 



2V A . k 2V Ak 



Z, 



C.k 



w A , k (z A M<i 2 + 4)] 



1/2 



± 



Z A , k tci 2 



(z 2 k + w A , k 5 k y/' 



+ o(|q| 4 ), (si) 



where 



and 



C,k 



± 



1 



2VA, k 2V A ,, 



(Z^ k + AV A>k 8 k ) 



1/2 



(52) 



When Zc,k 7^ 0, both modes have a quadratic dis- 
persion for small q. The modes w+(q) and w_(q) have 
positive and negative effective mass, respectively. Thus 
a; + (q) (w_(q)) corresponds to a particle-like (hole- like) 
excitation. At the transition to the superfluid phase 
(4=o 0), A fc =o,+ (Afc =0l -) vanishes if Z c , k =o > 
(Zc,k=o < 0), but the particle-hole excitation gap 
Afe =0 — Afc = o + — Afc =0 ,- remains finite. The critical 
mode energy u ~ q 2 being quadratic in q, the dynamical 
critical exponent takes the value z — 2. 

When Zc.k — 0, the excitation spectrum takes the 
particle-hole symmetric form 



w±(q) = ± 



/ z A . k tq 2 



1/2 



± w 



A 



where 



A fc = 



( 4 



Z A ,kt 



1/2 



1/2 



(53) 



(54) 



At the transition (4=0 = 0), the particle-hole excitation 
gap 2Afe = o vanishes and the dispersion w±(q) = ±Cfc=o|q| 
becomes linear, which implies that the critical dynamical 
exponent takes the value z = 1. We can also understand 
this result as a direct consequence of the (relativistic) 



(2) 

k—0 7 ij 



[Eq. (41)] when 
= 
1)- 



Lorentz invariance of the vertex T 
Zc k=0 vanishes. The quantum critical point 4=o 
then coincides with the critical point of the (d 4 
dimensional XY model. 

We conclude that the universality class of the 
superfluid-Mott-insulator transition depends on whether 
Zc,k=o{ n — 0) vanishes or not. ' The same conclusion 
can be reached from the superfluid phase by considering 
the spectrum in the limit no fc=o 0. 



E. Flow equations 

Since we do not have an explicit (approximate) form 
of the effective action T k [</>* , <j>] , we cannot directly use 
Eq. (7) to obtain the RG equations satisfied by 14 (n) 
and r a>k (n) (a = A,B,C). We can nevertheless ob- 
tain RG equations for the effective potential V k (n) and 

( 2) 

the two-point vertex 1^ (q;n) in a constant field within 
a simplified Blaizot-Mendez-Galain-Wschebor (BMW) 
scheme. 57,58 

Equation (7) leads to the RG equations 

diV k (n) = - i jf dtR k (q) [G M (9; n) + G fe , tt (<?; n)} (55) 



~ \ E 9iG k ^ ll2 (q;(l))r[ 4 l i2ii {p,~p,q,-q-(f)) 

~\ E { r S 2 i 3 ^'9 I -p-9; < / , ) r i 3 k ll (-^p+9,-9;0) 

x [diG k .i ll2 {q;(j))]G k ^ li {p + q;4>) + (jp O -p,i t+ j)\ 

(56) 

(/ = ln(fc/A)) for the effective potential and the two-point 



vertex in a constant field (p. G k — —(TV + R k )~ 
single-particle propagator. We use the notation 



1 



d d q 
(27r) d 



00 

-™ 27T ' 



is the 



(57) 



where the momentum integral is restricted to the Bril- 
louin zone ] — n, n] d . The operator di = {diR k )dn k acts 
only on the I dependence of the cutoff function R k . The 
BMW approximation is based on the following two ob- 
servations, i) For a given momentum q, the frequency 
integral in (56) is dominated by the region \u\ < wjjT(q) 
where wjjT(q) is the lowest excitation energy defined by 
the propagator G k . Since the function diGij(q; <j>) is pro- 
portional to diR k (q), the integral over the loop momen- 
tum q in (56) is dominated by values of |q| of the order or 
smaller than k. It follows that the important frequency 
range for the loop integral is \uj\ < uj k where uj^ is the 
typical value of oj^ (q) for |q| ~ k. In the superfluid phase 
10 k ~ c kk (cfe is the velocity of the Goldstone mode), 
while in the Mott insulating phase u> k can be deduced 
from (51,53). ii) Because of the cutoff function R k (q), 
the vertices r£ (<fr • • • q n ) are smooth functions of mo- 
menta and frequencies in the range |qj|/fc, \oJi\/oj k <C 1. 
These two properties allow us to expand the vertices in 
the rhs of (56) in powers of q 2 //c 2 and ui 2 /u^ 2 . To lead- 
ing order, one simply sets q = in the three- and four- 
point vertices in Eq. (56). We can then obtain a closed 

(2) . "7 

equation for T k by noting that 



rg>,-p,o;0) 



1 



Wd<l>i 
1 d 2 



9 rg^fe. 



42) 



(58) 



Furthermore, properties (i) and (ii) allow us to use 

(2) 

the derivative expansion of the two-point vertex T k 
[Eq. (41)] to obtain the propagator G k to be used in 
the RG equations (55,56). Since it is however crucial to 
retain the full lattice structure at the beginning of the 
RG flow (fc ~ A), we take 

T A . k (q;n) = Z A>fe (n)e q + V A , k {n)u 2 + V£(n), (59) 
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which coincides with (41) for |q| <C A. We have intro- 
duced the shifted dispersion e q = i q + 2dt (e q ~ tq 2 for 
|q| <C A). Following Ref. 3, we define Za ^(ti) as 



1 d 
ZA,k(n) = - lim — -KT A!k (q;n) 
t q^O c*q z 



(60) 



so that ZA,k{nQ^k) has the meaning of a field renormal- 
ization factor (and should not be confused with a renor- 
malization of the hopping amplitude between nearest- 
neighbor sites 3 ). For k <C A, Eqs. (59) and (60) are equiv- 
alent to the so-called LPA' approximation (LPA stands 
for local potential approximation). For k ~ A, these 
equations can be justified by noting that in this limit 
ZA,k{n) — ZA,K{n) = 1 so that approximating the renor- 
malized dispersion by -^A,/c( rl )£qj which is valid for small 
q when Za,u is defined by (60), is expected to remain 
approximately valid in the whole Brillouin zone. The 

(2) 

derivative expansion of the local vertex I\ c is further 
discussed in Appendix C. 

Although the final flow equations rely on a derivative 
expansion of the vertices, they cannot be derived directly 
from a simple Ansatz of the effective action T^. The rea- 
son is that it is not possible to propose an Ansatz for 
the initial effective action Fa[0*,0] [Eq. (16)] based on 
a derivative expansion, since we do not know its expres- 
sion for arbitrary time-dependent fields. In the BMW 
approach, we deal only with quantities computed in a 
constant field which, for k = A, can easily be obtained 
from the Hamiltonian in the local limit as explained in 
Sec. II B. 

As far as the momentum dependence of the vertices is 
concerned, our BMW approximation (supplemented with 

(2) 

a derivative expansion of T y k ' to obtain Gk) is as legit- 
imate as the original one. It is however more ques- 
tionable regarding the oj dependence. Contrary to the 
momentum integral, the frequency integral in Eq. (56) 
is not exponentially cut off by the regulator i?j,(q). 
The integrand typically decays as a power of l/|w| for 
\ui\ w_(q), so that the contribution of large frequen- 
cies is small but not negligible. The reason why the 
BMW approximation nevertheless leads to accurate re- 
sults (see Sec. Ill) can be understood as follows. In the 
weak-coupling limit, the frequency dependence of the ver- 
tices and rl 4 ^ is weak,''" so that setting the loop 
frequency uj to zero in T^, 3 ' and r£ 4 ' , as well as using 

(2) 

a derivative expansion of Ti should be justified. In 

the strong-coupling limit, r^ 3 -* and do depend on 
frequency but this dependence is controlled by U as in 
the local limit. Since U w_(q), except deep in the 
Mott phase where the strong-coupling RPA is already 
a good approximation to the k = results, it appears 
again justified to set the loop frequency to zero in three- 
and four-point vertices and use a derivative expansion to 
obtain the propagator Gk- The use of a cutoff function 
Rk(l) acting both on momentum and frequency would 
put the BMW approximation on a firmer basis, but 



such a cutoff function would be incompatible with the 
initial condition S + AS\ = 5i oc of the lattice NPRG. 61 
The numerical solution of the flow equations can be 
further simplified by approximating Va, k{n) and ZA,k{n) 
by V A ,k = VA,k(no,k) and Z A .k = ZA,k(no,k)- 6 To obtain 
an accurate description of the phase diagram, it is never- 
theless necessary to keep the full n-dependence of Zc k (n) 
and V/c(n). 10 The n-dependence of Zc,k(n) is also neces- 
sary for a good description of the critical behavior at the 
multicritical points (Sec. IV). Away from the multicrit- 
ical points, and when accuracy is not the primary goal, 
it is possible to approximate Zc,k{n) by Zc t k(no,k), and 
expand the effective potential to quadratic order about 
its minimum, 



V .k + ¥(n - n .k) 2 if n , k > 0, 



VkM = , 

1 Vn k + 6 k n+^n 2 



if n ,fc = 0, 



(61) 



where 5k and are defined in (42). The BMW equa- 
tions and their various approximations are detailed in 
Appendix D. 



III. PHASE DIAGRAM 

To alleviate the notations, we drop the subscript k 
whenever we refer to a k = quantity (e.g. uq = riQ^k=o)- 

For given values of t, U and u, the ground state can 
be deduced from the values of the condensate density Uq 
(?io > in the superfluid phase). To obtain thermody- 
namic quantities, it is sufficient to integrate the RG flow 
down to k - 10~ 5 . In Ref. 19, we have shown that by in- 
creasing the functional character of the NPRG equations 
(e.g. by retaining the full n-dependence of 14 (n) rather 
than using the truncation (61)), we observe a nice conver- 
gence of our results, which we therefore expect to be close 
to the exact ones. The most accurate results, obtained by 
keeping the full n dependence of Vk{n) and Zc,k(n) are 
shown in Figs. 5 and 6. Both in three and two dimensions, 
the transition line between the superfluid phase and the 
Mott insulator is very close to the QMC result: 9,10 the 
tip of the Mott lobe (t/U = 0.0339, (i/U = 0.3992) dif- 
fers from the QMC data only by (0.001%, 3%) in three 
dimensions, while in two dimensions the tip is located 
at (t/U = 0.060, [ijU = 0.387), which corresponds to 
a relative error of order (1.5%, 4%). For comparison, 
in Figs. 5 and 6 we also show the mean-field (or strong- 
coupling RPA) phase diagram as well as the one obtained 
from Dynamical Mean-Field Theory (DMFT). 64,65 

An important characteristic of the Mott insulating 
phases is the vanishing compressibility k = dh/d/i = 0. 
The expression (44) enables us to determine the boson 
density n k directly from n ,fc, Va,&, Z c .k and X k by in- 
tegrating dfik/dfi. The unknown integration constant is 
easily fixed since we know that the density vanishes for 
u = —2dt. Alternatively, one can use the fact that the 
density is integer in the Mott insulator. This method not 
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t/U p/U 



FIG. 5. (Color online) Phase diagram of the 3D Bose- 
Hubbard model. Only the first Mott lobe (n = 1) is shown. 
The (green) dashed line shows the mean-field (or strong- 
coupling RPA) phase diagram. The QMC data are obtained 
from Ref. !) and the DMFT data from Ref. 64. 




FIG. 6. (Color online) Phase diagram of the 2D Bose- 
Hubbard model. The QMC data are obtained from Ref. 10. 



only avoids to numerically compute dV^^/d^i (which re- 
quires to solve the RG equations for nearby values of /i) 
but also turns out to be much less sensitive to numerical 
noise. Figure 7 shows the density n as a function of the 
chemical potential /j, for various values oit/U and d = 2. 
The vanishing compressibility k = in the Mott insu- 
lating phase n = 1 is clearly visible. In the figure, the 
density is obtained from k and the condition h = 1 in the 
Mott phase. If we use the condition = —2dt) = 0, 
we obtain n = 1 ± 0.03 (n = 1 ± 0.045) in the three- 
dimensional (two-dimensional) Mott phase n = 1. The 
error is more pronounced near the tip of the Mott lobe. 



IV. CRITICAL BEHAVIOR 

In Sec. II D, we have seen that the universality class 
of the superfluid-Mott-insulator transition depends on 
whether Zc(n = 0) vanishes or not. We must therefore 
determine the value of Zc(n = 0) at the transition. 

Let us consider the parameter 8(t,fj,) = V'(0) as a 
function of t and /i (with U fixed). In the Mott insulator, 



FIG. 7. (Color online) Density n 
values of t/U and d = 2. 



the Ward identity (50) becomes 



rik=o vs fi/U for various 



Z c = Z c (0) = 



(62) 



Since the transition line is defined by 8(t, (j,) = 0, Zq 
vanishes at the tip of the Mott lobe (fi — fi c ) where the 
tangent to the transition line is vertical. Moreover, since 
8(t, n) > in the Mott phase, on the transition line Zq is 
positive (negative) for fx > fx c (< fx c ). We deduce that at 
the tip of the Mott lobe, the quantum critical point co- 
incides with the critical point of the (d+ l)-dimensional 
XY model. The dynamical critical exponent takes the 
value z = 1 in agreement with the Lorentz invariance of 
the effective action when Zq = (Sec. II D 2). The up- 
per and lower critical dimensions are therefore d+ = 3 
and d~ = 1, respectively. The transition from the Mott 
insulator to the superfluid phase is driven by the vanish- 
ing of the particle-hole excitation gap, while the density 
is conserved. ' The critical point is a multicritical point 
as two parameters (t/U and n/U) have to be fine tuned. 
Away from the Mott lobe tip, Zq is nonzero and the dy- 
namical critical exponent takes the value z = 2. This 
transition, which is driven by a density change, is mean- 
field like for d > 2 (with logarithmic corrections at the 
upper critical dimension e?+ = 2). 



A. Multicritical point 

The critical behavior at the tip of a Mott lobe can be 
understood from the linearized flow equations. If we set 
Zc,k(n) = 0, we recover the flow equations of the (d+ 1)- 
dimensional XY model with one relevant direction in the 
space of parameters of the effective action. The flow of 
the corresponding scaling field (which we denote by r) de- 
termines the exponent v. Since Zc ^(n) enters the propa- 
gators quadratically, it does not enter the linearized flow 
equations (except of course its own RG equation) . Thus 
Zc,k(no,k) corresponds to the second relevant direction 
of the flow and is orthogonal (in the parameter space of 
the action) to the critical surface. 
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multicritical point generic transition is however possible that the regulator i? fe (q), which does 

not satisfy the Lorentz invariance of the effective action 

q 1/^ at the multicritical point, modifies the expected critical 

q ( v A,k V^ 2 ^, I z c,k \ u behavior. Setting Zc,k(n) = in the flow equations, we 

\Z A ,k£k) \Z A ,kCkJ g n( j 

h k- d {V A , k Z A , k e k yl 2 n k~ d Z c , k n 

k ~ d {^Yy^) fc " d fehW W, k =VA, k -^ 2 . (66) 

S k (Z A ,k e k) S k (ZA,k e k) &k 

k d VXl /2 (Z A ,kek)~ 3/2 ^k ^(Zc.kZAMky 1 ^ Given the small value of the anomalous dimension in the 

Z (~) (V Z )~ 1/2 Z (n) (d+l)-dimensional XY model (d = 2, 3), the results r) A = 

c,k\n) ( A,k A,k£k) c,k{n) ^ z — i are nevertheless satisfied to a very good 

VA,k ZA,k^kZ^ 2 k VA, k accuracy (see below). 



TABLE II. Dimensionless variables (Zc,k = Zc,k(no, k ))- 



1 . 2D multicritical point 



The behavior of the system near the multicritical point 
(t c ,/z c ) is best understood by considering the singular 
part V s (r, Zc) of the effective potential (Zc = Zc(no)). 1 
For small r and Zc, and d < 3: 



V s (r,Z c ) = s- d -*V B (a 1 t v r,8Z c ) 
= \r\ v ^V s {\r\- v Z c ). 



(63) 



Here we anticipate that the eigenvalue related to the scal- 
ing field Zc is equal to one (see below). V s being finite 
and nonzero in the limits r — > and Zc — > 0, V s (x) 
must behave like a constant when x — > and like x d+z 
when x — > oo. Moreover, r and Zc are presumably an- 
alytic functions of t — t c and fi — fx c , and must vanish 
linearly with t — t c as we approach the multicritical point 
on a typical path (i.e. a path which is not vertical in 
the (t/U,/j,/U) plane). Since the critical exponent v of 
the XY model satisfies 1 — v > for all d + 1 > 3, the 
argument of V s in (63) vanishes as t — t c — > 0. Given 
that V s (x) — > const as x — > 0, we conclude that Zc drops 
out of the scaling relation (63) and the multicritical point 
looks like an ordinary XY critical point as shown explic- 
itly below by the NPRG results. 

To make the fixed point manifest when the system is 
critical, we use the dimensionless variables defined in ta- 
ble II. The anomalous dimensions are defined by 



r)A,k = -di lnZ Aik , 
Vv.k = ~di In V A ,k- 



(64) 



The dimensionless frequency variable cl> (Table II) allows 
us to define a (running) dynamical critical exponent Zk = 
[u] from \Z AX \ = -rj Atk , [V A . k ] = -r)v,k, and [&] = 0, 
which gives 



Zk = 1 - 



VA,k - Vv,k 



(65) 



Here [X] denotes the scaling dimension of the variable 
X (momenta having as usual scaling dimension 1). At 
the multicritical point, we expect rj A = rjy and z = 1. It 



Let us first discuss the two-dimensional case. We find 
that Zc(ri) vanishes for fi — 0.382, slightly away from 
the lobe tip located at fi = 0.387 (we now focus on the 
Mott insulating phase n — 1). We ascribe this slight dis- 
crepancy to the fact that the local gauge invariance (see 
Sec. IIC), which leads to the Ward identity (62), is not 
strictly satisfied in our approach since it is violated by 
both the BMW approximation and the derivative expan- 
sion. On the other hand, the fact that the multicritical 
point lies very close to the tip of the Mott lobe indicates 
that the local gauge invariance remains nearly satisfied, 
and all consequences discussed in Sees. IIC and II D ap- 
ply- 

Figure 8 shows the RG flow at the multicritical point. 
The plateaus observed for the dimensionless conden- 
sate density fiQ k and coupling constant A^, as well as 
for the (running) anomalous dimensions r] At k and T)v, k , 
are characteristic of critical behavior. We clearly see 
the emergence of the Lorentz invariance as k decreases: 
Zc,k(no,k) ~ k is suppressed while rjA,k and rjv.k become 
nearly equal (implying Zk — 1). We find the critical expo- 
nents v = 0.699, r\ A = 0.049, rj v = r) A (l - W 4 ) = 0-049 
and z — 1.000, to be compared with the best known esti- 
mates v — 0.671 and r\ = 0.038 for the three-dimensional 
XY model.'" The exponent v is deduced from the run- 
away flow from the critical surface when the system is 
nearly critical (e.g. tlq ^ — n^(x e~ l l v with the critical 
value of no). 

Figure 9 shows |Zc,fc("o,fc)| near the multicritical 
point. It first decreases towards zero as the multicriti- 
cal point is approached. Then \Zc,k(no,k)\ ~ 1/^ ~ e~ l 
increases as the flow runs away from the critical surface, 
so that the critical exponent associated with the scaling 
field Zc y k(nQ^k) is equal to one. The anomalous dimen- 
sions r] Aj k and T\v,k show the momentum range where 
the flow is controlled by the multicritical point, as indi- 
cated by the plateaus in t] A ^ and rjy^ in Fig. 9. The 
end of the plateaus determines the Josephson length 

O EE kj 1 = k~\ 

We now discuss the behavior of the system for /i = fi c 
and t — > tj. The condensate density must vanish with 
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FIG. 8. (Color online) (Top) Dimensionless condensate den- 
sity no,h, coupling constant Afe and Zc,k(no : k) vs ln(A/fc) at 
the multicritical point n = 1 for d = 2. (Bottom) Anomalous 
dimensions rjA.k and ryv,fc vs ln(A/fc). The inset shows that 
Eq. (66) is satisfied when k — > 0. 
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FIG. 9. (Color online) (Left) In Zc,k(no,k) vs ln(A/fe) near 
the multicritical critical point. (Right) Anomalous dimen- 
sions r\A,k and rjv,k- The end of the plateau determines the 
Josephson length £j = kj . 
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FIG. 10. (Color online) Condensate density no, superfluid 
stiffness p s , compressibility « and Mott gap A vs \t — t c \ near 
the multicritical point (t c ,p c ) [d = 2]. The crosses show the 
critical behavior (67-69) with v ~ 0.699 and r\ ~ 0.049. 
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FIG. 11. (Color online) Same as Fig. 8 but for a three- 
dimensional system. 
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FIG. 12. (Color online) Superfluid stiffness p s , compressibility 
k and Goldstone mode velocity c = J p a /tz vs t — t c (d = 3). 

the critical exponent 2/3 = z/(d + z — 2 + 77), 

no ~ (t-i c )" (d+ *- 2+,,) . (67) 

From the scaling dimension [p s ] = d + z — 2 of the su- 
perfluid stiffne ss an d the fact that the Goldstone mode 
velocity c = p s /n remains finite due to the Lorentz in- 
variance of the effective action in the limit k — > 0, we 
expect 



Ps ~ (t-t c y (d+z - 2 \ 

K^{t-t c )'^ d+Z - 2 \ 



(68) 



Equations (67,68) agree with the results obtained from 
the numerical solution of the RG equations (Fig. 10). 

In the Mott phase, since the gap has scaling dimension 
[A] = z, it must vanish as 



(t c -ty 



(69) 
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for t — > t c , again in agreement with the results obtained 
from the RG equations (Fig. 10). 69 



2. 3D multicritical point 

In three dimensions, the system is at the upper critical 
dimension (d 4- z = 4) and the transition is governed 
by the Gaussian fixed point with logarithmic corrections 
due to the marginally irrelevant coupling constant A&. 
The numerical solution of the flow equations show that 
X k ~ l/|lnfc|, Z c ,k{no,k) ~ fe/|lnfe| and f]A,k,riv,k ~ 
1/| In /c | 2 , while h ^ converges to its fixed point value Hq 
logarithmically (Fig. 11). Figure 12 shows p s , k and c as 
a function of t — t r . 



B. Generic transition 

For all transition points away from the lobe tip, the 
Va,/cW 2 term is irrelevant with respect to to the Zc,k^ 
term. The dynamical critical exponent is z — 2 and a 
simple dimensional analysis shows that the upper critical 
dimension is d+ = 2. The transition is therefore governed 
by the Gaussian fixed point (with logarithmic corrections 



for d = 2) defined by 
where 



A* 



V\ = and rf A = tc = 



Vc,k = -dilnZ c ,k(n ,k)- 



(70) 



The dimensionless variables used to study the generic 
transition are defined in Table II. Determining the dy- 
namical critical exponent as in Sec. IV A, we find 



Zk 



2 - Vam + Vc,k- 



(71) 



In three dimensions, linearization about the Gaussian 
fixed point gives 



4 - 
3n 0>fc + 3^2 Va 



din ,k : 

diVAk = 2VA,k> 



(72) 



and r)A,k — Vc,k — 0. We deduce that A^ ~ k, Va,r ~ k 2 
at the critical point, in agreement with the numerical 
solution of the flow equations (Fig. 13). Figure 13 also 
shows that riA.k,Vv,k ~ k 3 while the relevant variable 
no.fc vanishes linearly with k at the critical point. 

When a generic transition point (t G (/x),/i) is ap- 
proached on a path of constant chemical potential \i by 
varying t — t c (n), we observe the mean- field behavior 



p s ~t- t e (n), 
K ~ const, 



(73) 



for t —> t c {n) (Fig. 14). The compressibility n remains 
finite at the transition and the velocity c vanishes. 
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FIG. 13. (Color online) RG flow at the three-dimensional 
generic transition. The dotted lines show fits to no,fe ~ k, 
A fc ~ k, V A ,k ~ k 2 and T)A,h, Vc,k ~ fc 3 . 
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FIG. 14. (Color online) Superfluid stiffness p s , compressibility 
k and Goldstone mode velocity c vs t — t c (fJ,) near a generic 
transition point (t c (p),n) [p = 0.7U, t c (p) = 0.021917 and 
d = 3]. The dotted lines show the mean-field critical behav- 
ior (73). 
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FIG. 15. (Color online) Same as Fig. 14 but for the upper 
critical dimension dt = 2. The compressibility diverges when 
t -> U(p). 



At the upper critical dimension [d = d£ = 2), the 
mean-field behavior is corrected by logarithmic terms. 
The marginally irrelevant variable A& is suppressed as 
|lnfc| _1 , while the relevant variable no,fc vanish as 
|lnfc| _1 at the critical point. We observe a divergence 
of the compressibility n as the phase transition is ap- 
proached (t — > t c (p)) (Fig. 15). 
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V. SUPERFLUID PHASE 

In Sec. VA, we show that our approach reproduces 
known results in the weakly-correlated (dilute) limit. In 
the following sections (Sees. VB-VD), we discuss the 
properties of the two-dimensional superfluid phase. The 
three-dimensional superfluid phase is briefly discussed in 
Sec. VE. 



which is nothing but the inverse healing length of the 
superfluid: kh — ^ ■ The flow is governed by the Gaus- 
sian fixed point uq = A = for k 3> kh, and is driven 
away from that fixed point when k <^ kh due to the fi- 
nite boson density. For k 3> kh, we can approximate the 
RG equations of Z^,k, Z c .k, Va,u and \ k by Eqs. (74,75) 
to leading order in (fc^/fc) 2 = \k h no.k h /^k- Moreover, to 
this order, the variation with k of the condensate density 
is determined by the equation 



The dilute limit 



<9fc(A fc n ,fc) = 



(78) 



At sufficiently low density, we expect the lattice to be 
irrelevant and the system to behave as a dilute superfluid 
gas of bosons with mass m = l/2t. The thermodynamics 
of a dilute Bose gas is controlled by the zero-temperature 
fixed point governing the quantum phase transition be- 
tween the vacuum (i.e. the Mott insulating phase with 
n — 0) and the superfluid phase. "' This quantum criti- 
cal point is nothing but a particular case of the generic 
quantum critical point discussed in Sec. IV B. 

Let us first consider the vacuum limit fik = ^o.fc = 
and /i = —2dt. The RG equations (D6) give 



so that 



7 (vac) 



J A.k 



7 (vac) 
J C,k 



1. 



V 



(vac) 



A,k 



0. 



(74) 



in agreement with the fact that the single-particle (nor- 
mal) propagator is not renormalized: G ni k{q) = (ico+fj,— 

tq) -1 . 56 The coupling constant A^ is given by 
8irta 

, (vac) _ 

(d = 2) 



A 



1 — T^ka 

4-n-t 



(d=3), 



(75) 



ln(f)+C 

for k <C A, where C is the Euler constant and a the 
"s-wave" scattering length for the Bose-Hubbard model, 

1 1 



a= l**t/U + A (d ^< 



(76) 



(recall we take the lattice spacing as the unit length) 
with A ~ 0.1264. Equations (75,76) are derived in Ap- 
pendix E. The dimensionless coupling constant A^. vac ^ — 

\u k d /e% (see table II) vanishes for k — > when d > 2 
and the quantum critical point governing the transition 
between the vacuum and the superfluid phase is Gaus- 
sian. The logarithmic vanishing of X^ ac ^ in two dimen- 
sions agrees with d = 2 being the upper critical dimension 
(Sec. IV B). 

In the dilute limit, the finite condensate density n 0j fe 
can be ignored as long as e k 3> 2Afcno,fc.'" This defines 
the characteristic momentum scale 



kh 



Ai 



~ n 0,k h 



1/2 



(77) 



n 0,k h 



Aaw 0: / 
At,. 



(79) 



(see Appendix E). This equation allows us to relate the 
chemical potential to the coupling constants at scale kh- 
In the low-density limit, 

V A (n) = V loc (0) - (m + 2dt)n + ^n 2 + 0(n 3 ), (80) 



where we have used Aa = V^'(0) and Eq. (31) with n\ oc 
0. From Eqs. (79,80), we deduce 



/i + 2dt = A A n 0iA ^ \k h n a ,k h 



(81) 



for n + 2dt> 0. 

Below the healing momentum scale kh, the finite con- 
densate density cannot be ignored. The Bogoliubov ap- 
proximation amounts to ignoring any further renormal- 
ization as k decreases from kh down to 0, that is approx- 
imating the k = effective action by its value at k — kh- 
The initial value Aa and the first part of the RG flow 
(k kh) takes care of T-matrix renormalization of the 
coupling constant A^ (which is usually included in the 
Bogoliubov theory). This is what allows us to express 
the final results in terms of the scattering length a rather 
than the bare interaction U. We expect the Bogoliubov 
approximation to be valid if the ratio between the mean 
interaction energy per particle and the typical kinetic en- 
ergy, 



7 



Afc h n 

tfj 2/d 



Afe 
t 



(82) 



is much smaller than unity. ' Note that in this equation 
we consider the coupling constant A/^ at scale kh- We 
can then define the dilute (or weak-coupling) limit by 
the conditions 7 <§; 1 and kh -C A. The latter inequal- 
ity ensures that the characteristic length scale associated 
with superfluid behavior is much larger than the lattice 
spacing (thus making the lattice irrelevant as far as the 
superfluid properties are concerned) . ' 

We are now in a position to reproduce the standard 
results in dilute Bose gases. In three dimensions, since 
a < 0.31 [Eq. (76)], the coupling constant A& ~ 8wta is 
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FIG. 16. Crossover line between the dilute limit (7 <C 1 and 
k h <C A) and the strongly-correlated limit. 



roughly constant for fc« A. Given that kh <C A in the 
dilute limit, we deduce 



kh ~ \/8nan, 
Afc h a 8nat, 



(83) 



/j, + 2dt ~ 8natn, 

while the sound mode velocity takes the value 

c ~ v / 2tX kh n 0M ~ 2tVAnan. (84) 

Since ZA,k h — 1, the superfluid stiffness is given by 

p s ~ 2tn Qtkh ~ 2in. (85) 

We have used n 0j fc h = n to leading order in the gas pa- 
rameter 7 ~ an 1 / 3 of the three-dimensional dilute Bose 
gas. In the limit U/t < 1, where a ~ {J/87ri (the T- 
matrix renormalization of the coupling constant is negli- 
gible), one finds 



(86) 
/i + 2dt ~ f/n, 

c ~ VMJfi. 

The domain of validity of the dilute limit is shown in 
Fig. 16. 

In two dimensions, the logarithmic vanishing of a1 v 
plays a crucial role. One finds 




Ann 
In kh.a\ 



1/2 



Ann 



Ant 



In \/na\ 
Ant 



1/2 



(87) 



|lnfc/ja| I In -\/na I 



and 



fi + 2dt 



Antn 
I ln-\/Ea| ' 



87m 



In \pTLQ 



1/2 



in agreement with the results obtained by Schick for a 
dilute two-dimensional Bose gas. ' Note in particu- 
lar that the small parameter 7 ~ \/\\a.\f%a\. When 



-4irt/U j 



is exponen- 



U/t -C 1, the scattering length a ~ e 
tially small [Eq. (76)] and one recovers Eqs. (86). The 
dilute limit is then simply defined by kh ~ \J (U /t)n <C A 
(Fig. 16). 

Even in the dilute limit 7 « 1 and kh <C A, the Bo- 
goliubov theory breaks down at the Ginzburg scale k G 
(see the discussion in the Introduction). The latter can 
be estimated from the one-loop correction to the Bogoli- 
ubov approximation. Using the results of Refs. 14 and 41 
with the bare interaction U replaced by X kh to take into 
account fluctuations at momentum scales larger than kh, 
we obtain 



k G 



kh exp (— const/V na 3 ^j (d = 3), 
I In vna\ 



(89) 



For U < t, Eqs. (89) become 

k h exp \~const/^/fi(U/t) 3 / 2 ] (d = 3), 

"°~< U jk > (,= 2). < 90 > 

Thus, in the limit 7 <c 1, the Bogoliubov approximation 
remains valid in a large part of the momentum range < 
I qj -C kh where the spectrum is linear, and breaks down 
only when |q| -C k G <C kh- Thermodynamic quantities 
i n o,k, Ps,k, Cfc, etc.) are nevertheless insensitive to the 
Ginzburg scale and can be obtained from the Bogoliubov 
theory. 14 ' 41 



B. RG flows 

The RG flow is shown in Fig. 17 for t/U = 10 and 
n = 1 (d = 2). For these values of t/U and n, the initial 
condition is well approximated by the Bogoliubov theory, 



n ,A ~ n, p s . A ~2tn, c A a (2Utn) 
A A ~ U, Z c ,a 1, U a ,a 0, 



1/2 



(91) 



and the T-matrix renormalization of X k for k S> kh is 
negligible. The healing scale is deduced from the numer- 
ical results and 



Z A ,k h tk h - X kh n 0ikh , 



(92) 



which generalizes the definition (77) to cases where %A,k H 
may differ from unity. The thermodynamic quantities 
fk>,kj Ps,k and c k vary weakly with k and remain close to 
their Bogoliubov estimates (91). On the other hand the 
Ginzburg scale k G manifests itself by a strong variation 
with k of X k , Zc,k and VA,k- We determine k G from the 
inflection point in VVfe, 



d?V A , k \ 



h — k. c 



0. 



(93) 
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FIG. 17. (Color online) RG flow in the superfluid phase, 
t/U = 10 and n = 1 (d = 2). 



=J\ / — At/I/ 



/ — Cfc/CA 

0.5 W/ -/W,2« 



5 10 15 



r(va 
C,ifc 




5 10 

ln(A/fc) 



15 



FIG. 18. (Color online) Same as Fig. 17, but for t/U = 0.05 



ri .. \ >i < si . ..„ „£ \( vac ) ov ,^ 7 ( vac ) ;,. 
(n = no.fe = 0) is also shown. 



and n ~ 10 . The RG flow of A). and Z c k in the vacuum 
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FIG. 19. (Color online) Same as Fig. 17, but for t/U ~ 0.062 
and ra = 1. 



Both kh and fcg [Eqs. (92,93)] are in good agreement with 
k h ~ y/(U/t)n and k G ~ (l//t)As fc . 



In the (perturbative) Bogoliubov regime k G <C < A, 
Afc, ^c?,fc an( ^ ^A.fe remain nearly equal to their initial 
values (91). We therefore expect the k — propagators 



Gn(q) 
Gtt(q) 
Git(«) 



fa 



.£■2' 

q 

2C/n 



-^ 2 +£ 2 

a; 



(94) 



E 2 ' 
q 



for |q| 3> fcc which is the familiar Bogoliubov form with 
Eq = [e q (e q + 2[/n)] 1 / 2 the Bogoliubov excitation energy. 
The spectrum crosses over from a quadratic dispersion to 
a linear sound-like dispersion at the (healing) momentum 
scale kh. For k G <C |q| <C k h , £q ~ CA|q| with ca ~ 

(2f/ni) 1/2 - 

In the (nonperturbativc) Goldstonc regime k <C k G , 
Afe, ^c,fe ~ vanish with fc — > while Va,*; — Va,1c=o 
takes a finite value. This regime is dominated by phase 
fluctuations, and characterized by the vanishing of the 
anomalous self-energy T, an ^{q = 0) = \kno,k ~ k and the 
divergence of the longitudinal propagator (see Eq. (95) 
below). 32 ' 39 ' The k = propagators (15) are given by 



Gu(q) = - 
Gtt(q) = - 

G lt (?) 



1 



1 



2An 2n CV w2 + C V ' 

1 2n c 2 1 



Va(u 2 + c 2 q 2 ) ,o s a; 2 + c 2 q 2 



, (95) 



2\71 Va{u 2 + c 2 q 2 ) p s dp w 2 + c 2 q 2 ' 



for |q|/|w|/c <C k G , where we have used c = (Z^i/V^) 1 / 2 , 
k = 2h Va, p s = 2tZ A n and lim fc _y Z c ^ k /X k = dn /dfi 
(see Sec. II D). The longitudinal propagator is obtained 
from Gg t k(q = 0) = -l/(2A fc n ,fc) by replacing A fc ~ k 
with Ca/w 2 + c 2 q 2 . 14 In the Goldstone regime, the exis- 
tence of a linear spectrum at low energy is due to the 
(relativistic) Lorentz invariance of the effective action 
(Z G ,k while VA,k — > Va > 0) and not to the finite 
value of the anomalous self-energy S an (q — 0) as in the 
Bogoliubov regime. Quite remarkably however, the value 
of the sound-mode velocity is insensitive to the Ginzburg 
scale k G . These results agree with previous studies of 
interacting bosons in continuum models. 13 ' 14,17 ' 18 

In Fig. 18, we show the RG flow for t/U = 0.05 and 
n ~ 10~ 4 . Although t/U <C 1, the very small value of the 
density ensures that the system is in the dilute limit with 
k G < k h < A (l h ~ -3 and l G ~ -8). For k > k h , 
the flow of the coupling constant X k coincide with the 

flow in vacuum (Afc ~ A^. vac ^). For k G <C k <C kh, the 
variation of X k is weak: in the momentum range k G <C 
|q| "C kh, the behavior of the system is well described by 
the Bogoliubov theory but with renormalized parameters 
(Xk h , n o,k h , etc.). In the Goldstone regime k -C k G , we 
recover the infrared behavior discussed above. 

As t/U decreases (at fixed density n), the dimension- 
less coupling constant 7 increases and eventually becomes 
of order one. A typical flow in the strong-coupling regime 
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A* 



n h 



superfluid Z* A V% k k rig 

multicritical point k~ v fc"" k fc 1 " 2 " k 1+rl 

generic transition Z\ VX Z c |lnfc| _1 fc 2 1 In fe| 1 

insulator Z\ V% Z c A* 



r 



0.5 



— kU 

■ ■ - p s j'2tn 

n /n " 

c/V2tUn 

h/A - 
-■ ka/A 



TABLE III. Infrared behavior of the two-dimensional Bose- 
Hubbard model. The stared quantities indicate nonzero 
fixed-point values and n denotes the anomalous dimension 
at the three-dimensional XY critical point. Zc,k stands for 
Zc,k(no,k)- 



7 1 is shown in Fig. 19 for t/U ~ 0.062 and n = 1. 
There is no Bogoliubov regime any more, kh ~ kc ~ A, 
and the condensate density nn.fc and the superfluid stiff- 
ness p Sj fe, are strongly suppressed. 

It is instructive to compare the infrared behavior in the 
superfluid phase or the Mott insulating phase with the 
critical behavior at the superfluid-Mott-insulator transi- 
tion (Table III). Both in the superfluid phase and at the 
multicritical points, the infrared behavior is character- 
ized by a (relativistic) Lorentz invariance. 
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FIG. 20. (Color online) Condensate density no, superfluid 
stiffness p s , compressibility k, velocity c, and characteristic 
scales kh and kc vs t/U at fixed density n = 1 (d — 2). Black 
dots show fits kh oc \J nU jt and kc oc \Jn(JJ /t) 3 . 

(a) Weakly-correlated superfluid phase: y <SC 1 
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Bogoliubov 
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(b) Strongly-correlated superfluid phase: y » 1 



k G ~k h ~A 



Goldstone regime 



C. Characteristic momentum scales: kh and kc 

Figure 20 shows k G , k h , n , k, p s and c vs t/U at 
fixed density n = 1. We see a sharp crossover between a 
weakly-correlated (kc *C kh, ~ n and p s ~ 1th) and a 
strongly-correlated (kh ~ kc ~ A, no -C n and p s <§C 2tn) 
superfluid phase as t/£/ is decreased. Close to the mul- 
ticritical point, there is a critical regime where the flow 
of rin.fc, Afc, ?7A,fc and Tjv^ shows plateaus characteristic 
of critical behavior (Sec. IV). The critical behavior ends 
at the Josephson momentum scale kj (Sec. IV A), and 
for k <C kj we recover the Goldstone regime of the su- 
perfluid phase. The behavior in the superfluid phase at 
fixed commensurate density h (n integer) is summarized 
in Fig. 21. Except for lattice effects (which force kh, kc 
and kj to be at most of order A) we recover the behavior 
of the (d+ l)-dimensional O(N) model. 39 

Figure 22 shows kh and kc at fixed t/U for h varying 
between and 1. In the small density limit 7 < 1, our 
numerical results for kh and kc agree with Eqs. (77,87). 
For n ~ 0.5, kc and kh become of the same order and 
7 > 1 as expected for a strongly-correlated superfluid 
phase. The behavior near the Mott insulating phase 
(1 — n <C 1) is similar to the low-density limit and reflects 
the fact that the transitions from the superfluid phase to 
the vacuum or the Mott insulating phase h = 1 belong to 
the same universality class. In particular, near the tran- 
sition to the Mott insulating phase h = 1, the system is 
effectively in the weakly correlated limit (7 = A kh /t <C 1) 
in agreement with the fact that the quantum critical 
point is Gaussian for d > 2. The analogy between the 



(c) Critical regime 

kj k G ~ A 

I 1 1 *■ k 

Goldstone critical 
regime regime 

FIG. 21. Behavior of the superfluid phase vs momentum scale 
k at fixed commensurate density h. kh'- healing scale, fcc: 
Ginzburg scale, kj: Josephson scale. 



2.5 




0.2 0.4 0.6 0.8 
fl 



FIG. 22. (Color online) Condensate density no, superfluid 
stiffness p s and characteristic scales kh and kG vs density 
h for t/U = 0.04 (d = 2). Black dots show fits k h oc 
(ft/| In V / H) 1/2 and k G oc k h X kh /t [Eq. (89,87)]. 
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FIG. 23. (Color online) Normalized spectral weight Z of the 
sound mode [Eq. (98)] vs t/U (n = 1 and d = 2). 



FIG. 24. (Color online) RG flow in the three-dimensional 
superfluid phase for t/U — 0.034 and n = 1. 



limits n <C 1 and 1 — n <C 1 leads to interesting conse- 
quences which will be discussed elsewhere. 



D. Low-energy spectrum 

The knowledge of the infrared limit of the one- 
particle Green's function enables us to obtain the spectral 
function 14 

1 



A(q,w) = --ImG n (q, w + i0 H 



2tt 



Im [Gu(q, u + iO+) + G tt (q, u + i0+)] 

(96) 



in the low-energy limit. From Eqs. (95) we deduce 



A(q,w) 



n c 



[<y(w-c|q|)-*(w + c|q|)] 



2p s |q| 

sgn(o;) Q(u - c]q|) 
47rn C ^2 _ c 2 q 2 



(97) 



for |q|, |w|/c <C fee. In addition to the delta peak due 
to the Goldstone mode, the spectral function exhibits a 
continuum of excitations which is a direct consequence of 
the singularity of the longitudinal propagator G\\. While 
the sound mode extends up to |q| ~ k^, the continuum 
is observed only at momenta and energies |q|,|w|/c < 
k,Q. In the weak-coupling limit, where the lattice does 
not play an important role, these results are in complete 
agreement with Popov's hydrodynamic theory. > The 
latter gives C = (itn/crio) 1 / 2 so that the ratio of spectral 
weights carried by the continuum and the sound mode is 
extremely small in the weak-coupling limit. 14,80 

It should be noted that there is no qualitative differ- 
ence between the weakly- and strongly-correlated super- 
fluid phases regarding the low-energy single-particle spec- 
trum [Eq. (97)]. In the strong-coupling regime however, 
the continuum of excitations due to the singular longitu- 
dinal propagator is expected to extend up to momenta of 
order k(j ~ A (i.e. over most part of the Brillouin zone) 
and carry a significant fraction of spectral weight. This 



expectation is confirmed by the suppression of spectral 
weight of the sound mode as the ratio t/U is decreased 
at fixed density. Figure 23 shows the spectral weight 
-Zq = nQc/2p s \q\ of the sound mode for n = 1, normal- 
ized by its value in the weakly-correlated limit i > [/, 



2|q|2, 




(98) 



Z remains close to one in the weakly-correlated super- 
fluid phase but is strongly suppressed in the strongly- 
correlated regime. It vanishes at the transition to the 
Mott insulating phase (t — t c ) with a critical exponent 
2/3 - u{d + z - 2) = VT), 



(t - t c r 



(99) 



for t -> t+. 

Being equivalent to the strong-coupling RPA, the ini- 
tial effective action Ta predicts the existence of a gapped 
mode in addition to the sound mode. 48-50 ' 79 ' 81 Whether 
this gapped mode is a true characteristic of the spec- 
trum (which would then show up in the propagator 
Gfc = o(<7; ^o)) is an interesting question which however re- 
quires a more refined NPRG analysis 15 ' 14 ' 1 '' 18 beyond the 
derivative expansion [Eq. (41)]. 



E. The 3D superfluid phase 

In the dilute limit 7 <C 1 (Sec. VA), the initial condi- 
tions of the RG flow are given by 



no,A — n , Ps,A — 2tn, ca ~ 2iv 4nan, 
A A ~87rai, Z CiA (n) ~ 1, V a ,a(ti) ~ 0, 



(100) 



and reproduce the Bogoliubov approximation. The flow 
for k < A is logarithmic and therefore very slow. This 
explains why the Ginzburg scale is exponentially small 
in the dilute limit [Eq. (89)] and irrelevant for most pur- 
poses. In the infrared limit, Gtt and Git are given by (95), 
while the longitudinal propagator 



Gu(?)~ln 



ck h 



\Juj 2 + c 2 q 2 



(101) 
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diverges logarithmically. 

A typical RG flow in the strong-coupling limit, where 
khikG ~ A, is shown in Fig. 24 for n = 1. Thermody- 
namic quantities very rapidly converge to their k = 
values. On the other hand, the flow of Zc,k{no,k) ~ 
|lnfc| _1 and Va,A: is logarithmic. 



VI. SUMMARY AND CONCLUSION 

We have presented a detailed NPRG study of the Bose- 
Hubbard model. Although we have only considered the 
zero-temperature limit, it is straightforward to extend 
the analysis to finite temperatures. The lattice NPRG 
seems to be the only available technique which treats 
fluctuations at all length scales on equal footing: 

• the lattice NPRG takes into account on-site cor- 
relations which are responsible for the very exis- 
tence of the superfluid-Mott-insulator transition. 
It is exact in the local limit (vanishing hopping 
amplitude), the latter corresponding to the initial 
condition of the NPRG at the microscopic scale 
A = V2d. 

• the lattice NPRG also takes into account critical 
fluctuations at the superfluid-Mott-insulator tran- 
sition. In this respect it is very similar to the stan- 
dard implementation of the NPRG in continuum 
models, the cutoff function -Rfe(q) playing the role 
of an infrared regulator. 

• as already known from previous studies in contin- 
uum models, 12 ~ 15 ' 17,18 the NPRG is a method of 
choice to study the superfluid phase. It is free of 
infrared divergences, satisfies the Hugenholtz-Pines 
theorem, and is able to describe both the (perturba- 
tive) Bogoliubov regime kQ <C k < A and the (non- 
perturbative) Goldstone regime k <g; kc- The lat- 
ter is characterized by a vanishing anomalous self- 
energy and a diverging longitudinal propagator. In 
the strong-coupling limit where kQ ~ kf, ~ A, there 
is no Bogoliubov regime and the whole RG flow be- 
comes nonperturbative. 



Our results agree with known results on the Bose- 
Hubbard model. In particular we reproduce the phase 
diagram obtained from QMC calculations with a typi- 
cal accuracy of 1-3 % and at a very modest numerical 
cost. 8,5 Moreover, we recover the two universality classes 
of the superfluid-Mott-insulator transition. The lattice 
NPRG enables a detailed study of the critical behavior 
near multicritical or generic transition points, which con- 
firms the original predictions of Fisher et al. 1 based on 
scaling arguments. 




FIG. 25. (Color online) Ground state energy Eo and occu- 
pation number n vs fi/U in the local limit (with vanishing 
external source). 
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Appendix A: Effective potential in the local limit 

In this appendix, we discuss the solution of the local 
Hamiltonian (18) when the ground state is degenerate for 
vanishing external sources (n/U integer). This degener- 
acy has important consequences for the effective potential 
V loc (n). 

When J* = J = 0, the Hamiltonian is diagonal in the 
basis {|m)} (Sec. II B). The ground state is the vacuum 
state |0) for /i < 0, and \m) for m < fi/U < m + 1. 
The ground state energy Eq and occupation number n 
are shown in Fig. 25 as a function of n/U. There is a 
quantum phase transition whenever fi/U is integer due 
to a level crossing. When \ijU — m (m integer), the 
states |m) and |m + 1) are degenerate. 

For an infinitesimal external source and (i/U = m (m 
integer), it is sufficient to consider the degenerate states 
|m) and \m + 1) to determine the ground state of the 
Hamiltonian (18). In this subspace, 



H 



m —J*\/m + 1 

-J\Jm + 1 e m+ i 



(Al) 



with e m = e m+ i - 2 
find the two states 

H = T2 

1 



u m(m + l). Diagonalizing (Al), we 



(\m) +e w \m + l)) , 



(A2) 



H = ~{\m)-e ie \m + l)), 



with eigenvalues 



E- = e m - \J\y/m + 1, 
E+ = e m + \J\y/m + l, 



(A3) 
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FIG. 26. (Color online) Superfluid order parameter 4>(J) = 
(6) (top) and effective potential Vj oc (n) (bottom) in the local 
limit. The source J is taken real. 



where 9 denotes the phase of the complex source J = 
| J\e l6 . The occupation number in the ground state |— ) 
is 



(— ri — > = m - 



1 

2' 



(A4) 



while the superfluid order parameter 



-|S|-) = T Vm + l 



(A5) 



is finite. We conclude that the U(l) symmetry is spon- 
taneously broken whenever ji/U is integer, 



lim lb) ^ 



(A6) 



although (b) = for | J\ = 0. 

Figure 26 shows the superfluid order parameter 0(J) 
obtained from the numerical diagonalization of the full 
Hamiltonian (18). For n = 0, we find \<j)(J = 0+)| = 1/2 
in agreement with Eq. (A5). The effective potential 
V\ oc (n) takes the usual form in a system with a spon- 
taneous broken symmetry, with a flat part (V^ c (n) = 0) 
for n < 1/4. By contrast, for \ij~U = 0.02, the superfluid 
order parameter cf>(J = + ) vanishes and V{ oc (n) > for 
all values of n. 




FIG. 27. (Color online) Derivative V{ oc (n) of the local effective 
potential for various values of fj,. The dotted lines show the 
large-field limit (B5). 



Appendix B: Large-field limit of the local effective 
potential Vi oc (n) 

The large-field limit of the effective action ri oc [</>* , <j> 
is obtained by considering the partition function 

Z loc [J*,J} = [ V[ip*, ^] e -S^[r M+I e dr(r,p+Ji,) 



(HO 

for | J | — > oo. In this limit, we expect the field to weakly 
fluctuate about its saddle-point value ip c defined by 



6S, 



loc 



J(t). 



(B2) 



Let us compute the effective action by including Gaussian 
fluctuations about the saddle-point solution ip c (one-loop 
order). The calculation is standard and gives'^ 1 



Tloc 

where 
G7 1 



= S loc [cf>*,<l>} + l r ftlnG, 



(B3) 



t,t : 
( d T 



,4>] = -6(t-t') 
M + 2[/|0(t)| 2 



d T ,- t i + 2U\0(T)\ 2 



(B4) 



is the inverse classical (local) propagator. By performing 
the trace in (B3) for a time-independent field </>, 8j we 
easily obtain the effective potential 



V\ oc (n) 



U 2 

[in + —n 



+ l - { [(/i - 2Un) 2 - U 2 n 2 } 1/2 + y. - 2Un} 

(B5) 



fin + ^-n 2 + O{n ), 



where 



p, = n + U 1 - 



V3 N 



(B6) 
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FIG. 28. (Color online) [V loc , A (iu; n) - V{ oc (n)]/U vs u/U for 
various values of n. fi — 0.2U (left) and fi — (y/2 — l)U (right). 
The dotted lines show the derivative expansion Va(«-)cj 2 . 
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FIG. 29. (Color online) Fi OCt c{iu;n) /U vs u/U for various 
values of n. The dotted lines show the derivative expansion 
Z c (n)u. [x = -0.2(7, 0.2(7, (y/2 - l)U and 0.6(7 (from top 
left to bottom right). 



derivative expansions 



V A ,\oc(n)uj 2 + V{ oc (n), 
Zc,\ oc {n)uj. 



(CI) 



The derivative expansion is remarkably accurate when- 
ever the chemical potential is negative or the condensate 
density large. In both limits, T\ OCt A{iw;n) ~ ^o C ( n ) an d 
rioc,c(* w ; n) ~ uj. Since a negative /i or a large n cor- 
responds to a system deep in the supcrfluid p hase, we 



l- 
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. \ v. 

N\'-. 



0.5 - 



/ 



■ 



H = -Q.2U 
H = 0.2U 
H={s/2-l)U 
H = 0.6 U 







FIG. 30. (Color online) Va,1oc(^) an d Zc.ioc(n) vs n for vari- 
ous values of the chemical potential /i. 



To one-loop order, the effective potential is given by the 
microscopic action S\ oc with a shift of the chemical poten- 
tial. It is straightforward to verify that higher-order con- 
tributions (e.g. those coming from two-loop diagrams) 
are at most of order O(n ) in the large-field limit. Equa- 
tion (B5) is in very good agreement with the numerical 
calculation of V\ oc (n) (Fig. 27). 



conclude that the derivative expansion is fully justified 
in this limit. 

More generally, we see that the derivative expansion is 
always valid in the limit \uj\ <C U. As argued in Sec. HE, 
except deep in the Mott phase (where the strong-coupling 
RPA is a good approximation to the k — results), U 
is a very large energy scale in the strong-coupling limit, 
and the knowledge of the vertices at energies |w| -C U 
is sufficient to solve the flow equations. We therefore 
expect the derivative expansion to be justified also in the 
strong-coupling limit. 

Figure 30 shows Va,1oc(^) and Zc.ioc(n) for various val- 
ues of the chemical potential \i. 



Appendix C: Derivative expansion of the local 



vertex Y 



(2) 



Figures 28 and 29 show the local vertices Fi OC j 4 and 
TiocC (also shown in Figs. 3 and 4) together with their 



Appendix D: Flow equations 

The flow equations in the BMW scheme can be found 
in Appendix C.l of Ref. 14. When ZA,k{n) and VA,k(n) 
are approximated by their values at the minimum no.k of 
the effective potential, the equations simplify into 
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diV k (n) = - - I diR k (q)[Gk,n{q;n) + G k ,tt(q;n)] 



(Dl) 



diZ C: k(n) = - 2nV k (n) 2 d U! [3 Jn,it(g, n) - 3J\ t ,n(q,n) + Ji ttit (q,n) - J tt ,it (q, «)L =0 



4nV fc (3) (n)14"(n)9 u [Jn,i t (q,n) - J lt>11 (g, n)] q=0 - An 2 V k ^ {n)Z' Cjk {n) Ji yl (g, n) 

^ fe (n) [Z tt (n) + J u (n)] - h{n)nZ'^ k {n) - 2nV k "(n)Z' c>k (n)J im (0,n) 
6nVl'{n)Z' c ^(n)J n ,n(0,n) 

n Q,k>^dl [Jn : tt(g) + 2Jit,it(?) + Jtt,\\(q)] q=0 + 2n Qtk Jn,\\(0)Z' C j t (n , k ) 2 
4n a ^ k X k Z / C j 1 (n 0tk )d UJ [J\\,it(q) - ^lt,ll(«)] q=0 , 



2 T /(3)/ 



7]A, k 



2\ 2 k n ^ k Z A , k e k k 2 -j + ^ 6»(e fc - e q ) 



1 1 
D 2+ D 2 I 



where 



D- = (Z A . k e k + V A ^ k u A ){Z A , k e k + V Aik u 2 + 2n . k X k ) + Z C M(n . k ) 2 uj 2 , 
D+ = [Z A>k {Adt - e k ) + V Ak LJ 2 ][Z A , k (4dt - e k ) + V A;k 0J 2 + 2n , k \ k ] + Z c , k ( n o,k) 2 



(D2) 

(D3) 
(D4) 

(D5) 



Z Aik = Z A>k (n 0>k ), V A , k = V Ajk (n ,k), Vam = -kd k \nZ Atk , and I = ln(fc/A). 

When Zc, k {n) is approximated by Zc }k {jiQ }k ) = Zc }k and 14(n) truncated to quadratic order [Eq. (61)], we obtain 



diV , k 



d,i? fc (q) [G k ,n(<l) + GkM] 



3 - 1 - 



<9;<5fc 



-2Afc7 n if n ,fc = 0, 
-A|[9J u ,n(0)-6J lt)lt (0) + J tt ,tt(0)], 
d 



diZ c . k = 2\ln 0jk —[Jtt,it(q) - JitMq) ~ 3 ^n,it(<?) + 3Jit,n(<?] q=0 



diV A . k = -2\ln ,kg^[Ja,tt(q) + Jtt,u(l) + 2 ^iut(g)] 9=0 , 



(D6) 



I 

with rj A>k given by (D4). We have introduced the coefR- effective potential, this gives 
cients 



I a (n) = / diG k , a (q;n), 



J al 3(q;n) 



[diG k . a (q';n)]G ky p(q + q';n), 



(D7) 



defined as — (T k + R k )~ 1 , with T k approximated by its 
derivative expansion (41,59). With Z A k (n) and V A}k (n) 
approximated by their values at the minimum no. k of the 



T A<k {q; n) = Z A . k eq + V A . k uj 2 + V fe '(n), 
r Btk (q;n) = V k "(n), 
^C,k{q;n) = Z c ,k{n)ui 



and 



where a, f3 = 11, tt, It. To alleviate the notations, we have 
omitted the subscript k in I a and J a p. The notation 
la, Ja/3 and G means that these quantities are evaluated 
for n = no.fe. The Green functions G k . a in (D7) are where 



G k .\\{q;n) 

Gk,u(q;n) 

G ki \t{q;n) 



(D8) 



F A:k (g;n) +R k (g) 
D k (q;n) 

T Atk (q; n) + 2nT Btk (q; n) + R k (q) 



D k (q;n) 



Tc,fc(g; n) 
D k (q;n) 



(D9) 



D k (q;n) = [T A:k {q;n) + R k (q)} 2 + 2nT Btk (q;n) 

x [T A>k (q; n) + R k (q)] + V c , k (q; n) 2 . (D10) 
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Lattice regulator 



Introducing the lattice density of states 



The lattice cutoff function (4) differs from cutoff func- 
tions used in the continuum, in particular due a sym- 
metric treatment of the low- and high-energy parts of 
the spectrum (Fig. 1). We can rewrite -Rfc(q) in the form 



where 



and 



•Rfe(q) = -^A,feefeSgn(f q )y q r(y q 



r( !/ ) = L_Ee(i_„) 



if tq < 0, 

■'let = \ tdt - e q , 

^ if t q > 0. 



(Dll) 



(D12) 



(D13) 



This gives 



^-Rfc(q) = Z A ^ k e k sgn{t c ^y CL [r) Atk r(y CL ) + 2y q r'(y q )] 

= Z j4 ,A : e fc sgn(t q )e(l - y q )[?M,fc(l ~ J/q) ~ 2]. 

(D14) 

i?fc(q) enters the flow equations always in the combina- 
tion 

r A , fe (g; n) + R k (q) = Z A , fc e q + V A , k uP + V k \n) + R k (q), 

(D15) 

where 

Z A , fe e q + E fe (q) = { ^' fe£fc * £ ^ £fe ' (D16) 
I ZA,fc£ q if e q > £fc 



for tq < 0, and 



^ D / x I Z A<k (Adt- ejfe) if e q >4di-e fe , 

[ ^4,fc e q it e q < 4dt - e fc 

(D17) 

for t q > 0. Equations (D14,D16,D17) lead to a significant 
simplification of the coefficients I a (n) and J a p{q\n). A 
typical contribution to / a or J a p reads 

F k = [ /'^flfc(q)/(Z A)fc e q + J2 fc (q),«) ) (D18) 

where / is a product of propagators G k — — (T k 2 '+R k )~ l . 
Since diR k (q) restricts the momentum integral to the 
domain y q < 1 where Z A ,k£q + -Rfc(q) is independent of 

q, 

F k = [ Q{-tq)diRM) I f{Z AM e k ,uj) 

J q J uj 



+ / 9(t q )fti? fc (q) / /(Z Aifc (4dt-e fc ).w). 



P(e)= / <S(e-e q ) (0 < e < Adt), 



(D20) 



we obtain 



/ e(-t q )9ii? fc (q) = / deZ>(e)0|fl fc (q) 

Jq JO 



-Z A ,fcefe / deT>(e) 
Jo 



VA,k 1 



(D21) 



and 



6(t q )5 z i? fe (q) = / deV{e)diR k {q) 

q J2dt 



4ri( 



J 4dt— e k 



VA,k 1 



Adt-e 



2 . 

(D22) 



Since the hypercubic lattice density of states is symmet- 
ric, T>(e) — V{Adt — e), the last equation gives 

f e(t q )9 i ii fc (q) = - f e(-t q )a ; i? fc (q) (D23) 

which enables us to rewrite (D19) as 
F k = f e(-t q )5 ; i? fc (q) 

Jq 

x / [/(^ ifc e fc ,w)-/(Z A)fc (4dt-e fc ),w)], (D24) 

with the momentum integral given by (D21). 

In the limit k <C A (or e k <C 2dt), the function 
/ \Z A,k(Adt— e k ) , uj) involves propagators with a large gap 
and is therefore negligible with respect to f{Z A>k e k ,uS); 
the flow is then governed only by the low-energy modes. 
Moreover, for e q < e k <C 2dt, the lattice does not matter 
and we can approximate e q ~ tq 2 , which leads to 

V{e) = 2v d "-^- (e«2tft) (D25) 



and 



e(-t q )Qi2 fc (q) = 8-j Z M k d 1 

a \ a + 2 



(D19) 



(D26) 

,/-■ - V ,,-■!) J 

(D27) 

which is the usual form for models in the continuum limit 
with the theta cutoff function. 



where v^ 1 = 2 d+1 n d ^T{d/2). We then obtain 
F k = % V -±Z A , k e k k d (l VA ' k 
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Appendix E: The vacuum limit 



and 



1. Scattering length on the lattice 

Let us first recall how the s-wave scattering length a 
is computed from the low-energy behavior of the T ma- 
trix in the continuum. For a contact interaction U, the 
retarded T matrix is defined by 



T*(w) U 



(El) 



where 



q 2e q - u - W 



q 2e q - u 



+ in / 5(w-2e q ) (E2) 



(V denotes the principal part) and e q = q 2 /2m is the 
dispersion of the free bosons. In three dimensions, 



T«(q 2 /m) = ^-^-_- (| q |_> ), (E3) 



m 1 + i|q|o 
while in two dimensions 



T«(q 2 /m) = - 



27r/m 



In 



|q|a 



C-i§ 



(|q|->0), (E4) 



where C is Euler's constant. Equations (E3) and (E4) 
define the s-wave scattering length a in three and two 
dimensions, respectively. 

A scattering length can be defined similarly in the 
Bose-Hubbard model. At low energy (e q <C t) the lattice 
does not matter and one can approximate the boson dis- 
persion e q = tq + 2dt by tq 2 . The bosons then behaves 
as free particles with an effective mass m = l/2f. Thus, 
equations (E1-E4) allow us to define a scattering length 
provided that we replace m by 1/24. 



a. d = 3 
In the low-energy limit < ui <C t, 

H R (u;)~ I ^~+ ™ / <5(cj-2e q ) 

Jq ^ e q Jq 

23(e) 7T fUJ\ , . 

j +^y, ( E5 ) 



q ^q 

8t 

o 2e 



where 23(e) = J q <5(e — e q ) is the density of states of the 
cubic lattice. The last integral in (E5) can be computed 
numerically while 23(u;) ~ v / a;/47r 2 t 3 ' /2 for small us. This 
gives 



1 A 



,1/2 



T r {uj) U t 



(E6) 



1 + ilqla 



for q — > with 



T R (2t\q\ 2 ) 8-n-ta 



1 1 

8nt/U + A 



(E7) 



(E8) 



and A ~ 0.1264. 



b. d = 2 

To compute H R (ui) in two dimensions, we use 
1 „ f ( 1 1 



q 2e q - a; 



= 7> 



2e q — u> 2iq 2 — CJ 
1 

2iq 2 - us ' 



(E9) 



Since the first integral in the rhs of Eq. (E9) is convergent, 
we can set u> — 0, which gives 



1 1 
2e^ ~ 2tc? 



1 W, 8 

— -2 G + - In -2 
2tir 2 V 4 7r 2 



(E10) 



where G ~ 0.916 is the Catalan constant. As for the last 
integral in Eq. (E9) we obtain 



q 

1 , /2*7T 2 

= In 

8nt \ UJ 

for < lu <C t. Since 



^flu2- 2 G 
4?rf 



i7r / 5(u> - 2e q ) ~ in I 5(uj - 2iq 2 ) 



we finally obtain 
1 



T R (2t\ q \ 2 ) 
for q — > with 



1 

~4~7rf 



2 

8t' 



\nl + C -i 



2V2 



-i-Kt/U-C 



(Ell) 



(E12) 



(E13) 



(E14) 



2. Coupling constant Xk 

The low-energy limit A: <C A of the coupling constant 
in vacuum (n/- = no,fc = 0) can be expressed in terms 
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of the scattering length a. X k can be obtained from the 
RG equation (D6) (see Sec. E 3) or more simply from 



1 1 1 



1 



A fc U 2 7 q e q + i? fc (q)' 



(E15) 



where e q + -Rfc(q) is given by Eqs. (D16,D17) with Za.u — 
1. For k -C A, we can ignore the effect of the cutoff 
function on the high-energy part of the spectrum (i q > 
0), 



1 1 f e(e k - e q ) , 1 f 6(e q - e k ) 



X k U 2 



c k 



2 / q e q 



(E16) 



Using 



a. d = 3 



1 /" e(e fe - e q ) fc 



2 ./q £fc 127T 2 t 



(E17) 



and 



1 /■ 9(e q - e k ) = 1 I 1_ _ 1 f 6(e fc - e q ; 
2i q e q 2i q e q 2 i q e q 

1 /" 1 fc 



2 ,/ q e q 47T 2 t 



(E18) 



(we have approximated e q ~ tq 2 for e q < £t <C A), we 
finally obtain 



1 1 1 f 1 fc 



A fc tf ' 2./ q e q 6^ 2 i 
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where a is the three-dimensional scattering length (E8). 
6. d = 2 

In two dimensions, we rewrite Eq. (E15) as 
111/" 1 



Afc C/ 2 7 q iq 2 + i? fe (q) 

1 /" ( 1 1 



2 7 q Veq + ^fe(q) tq 2 +i? fc (q) 



(E20) 



where i?fc(q) is obtained from i?fc(q) by replacing e q with 
tq 2 . Since the last integral in (E20) is convergent for 
k — >• 0, we can set fc = and use the result (E10). The 
first integral in (E20) can be expressed as 

1 f 9(fc-|q|) If e(7T- |q|)e(|q| -fc) 
e k 2 7 q tq 2 

1 /■ e(|q| - tt) 

2 y q 

= _^ + J-ln^ + J- (l n 2-V). (E21) 



We deduce 
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(E22) 



where a is the two-dimensional scattering length (E14). 

3. RG equation dkXk 

In the vacuum, Z A ,k = Zc.k = 1, Va,/c = and n ,fc 
0. the two-point vertex is defined by 

^A,k(q) = e q , ^BM(q) = Afc, r c ,fc(g) = w. (E23) 
The RG equation satisfied by Afc takes the simple form 



d k X k 



where D = [e c 



-A? 



2 / dkRk( q) j 
fc / £,3 



{22 W 2 [e q + i? fc (q)] 



-lOleq + iT-^q)] 3 }, 
i? fe (q)] 2 + cj 2 . We deduce 
A 2 /• a feJ R fe (q) 



2 J q K + ^R^q)] 1 



(E24) 
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and 
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q v e q + i? fe (q) e q + i? A (q) 



2 7 q e q + i?fe(q) 4dt 



(E26) 



Aa can be computed from the local action S\ — S\ oc in 
the vacuum (/U = —2dt), 

h = b + f u G M G w = b + **> (E27) 

where G(iuj) = (iu) + is the local (normal) propaga- 
tor in vacuum. From Eqs. (E26) and (E27), we recover 
Eq. (E15). 



4. RG equations in the dilute limit for k ^> kh 

We now consider the RG equations at finite density 
but in the dilute limit (kh <C A) for fc » fc^. To leading 
order in X k n 0jk , At satisfies the equation (E24). To ob- 
tain the equation satisfied by no tk , we must expand the 
propagator to first order in \ k no^ kl 



r< , . e q + -Rfc(q) . 2X k n ok , . ., 2 

Gn(g) = — - — p; + ^ [eq + J4(q)J , 



^ , s e q + fl k (q) 2X k n Q . k 2 
l^ttVl) — ~ ^ — w 



(E28) 
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where D is defined in Sec. E3. We will show below that 
the flow of Zam, Zc,k and Va,u leads to higher-order 
corrections. We can therefore set Za,Ui %c,k and Va,h to 
their vacuum values. This gives 



^^2A fc n , fe {5 W 2 [ eq + i? fe (q)] 
-3[e q + i4(q)] 3 }. (E29) 



Eq. (78) follows from (E24) and (E29). 

Let us now show that the flow of Za,Hi Zc,k and Vam 
give subleading contributions to Eq. (E29). To leading 
order in A^no^, the RG equation of Zc t k reads 



diZ c ,k = -A fc n , 



\2 



diRkjq) 
q L e(q)+ii fc (q)]3 

e(-t q )9 ; i? fe (q) 



1 



(E30) 

To estimate the order of magnitude of diZc.k, we can 
ignore the term (4dt — £fe) -3 (which is smaller than e^ 3 ) 
and use the approximate density of states (D25). This 
gives 



diZ 



l^C,k 



-s v f t \ k klk d -\ 



(E31) 



where we have used Afcno,fc — 
tegrating this equation between k 



tk\ to leading order. In- 



A and k <C A, we 



obtain 



At-fc 



k K h 



tk A ~ 



(E32) 



(ignoring the dependence of A& on k). In the dilute limit 
and for k > kh, the rhs is always small. For example, in 
three dimensions, one finds Zc,k h — 1 ~ k^a -C 1. 
Similarly, we obtain 



diV A ,k 



and therefore 



3 X2 



diR k (q) 



[e( q ) + Rk(q)Y 



Va„ 



X k k 



kK h 



t 2 k e - 



(E33) 



(E34) 



for k -C A. Since the frequency integrals extend up to 
uj oj tk 2 while e q ~ tk 2 , the Va,/cW 2 term in the propaga- 
tors can be neglected if VA,k <C l/tk 2 . Equation (E34) 
shows that this is indeed the case. For example, in three 
dimensions tkf^VA.kh ~ *C 1. 
Finally 



(E35) 



where we have calculated the integral over q in (D4) using 
e q = tq 2 . As Zc : k, ZA,k remains close to unity for k > 
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